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Abstract. We investigate some aspects of the ageing behavior observed in the 
contact process after a quench from its active phase to the critical point. In particular 
we discuss the scaling properties of the two-time response function and we calculate 
it and its universal ratio to the two-time correlation function up to first order in the 
field-theoretical e-expansion (e = 4 — d). The scaling form of the response function 
does not fit the prediction of the theory of local scale invariance. Our findings are in 
good qualitative agreement with recent numerical results. 
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The contact process (CP) has been introduced long ago as a simple model for the 
spreading of diseases [l] (see also El IH |H] ) . It is defined on a d-dimensional hypercubic 
lattice in which each site x is characterized by an occupation variable n x and is either 
active (ill, n x = 1) or inactive (healthy, n x = 0). At each time step of the dynamics 
a site x is chosen at random. If x is inactive then it becomes active with a rate given 
by A times the fraction of its neighbouring sites which are already active (infection). 
Otherwise, if originally active, x becomes inactive with rate 1 (healing). In the long-time 
(non-equilibrium) stationary state all lattice sites will be eventually inactive (absorbing 
phase) if A is small enough, whereas for large enough A a finite fraction of sites will be 
active (active phase). These two phases are separated by a continuous phase transition 
for A = A c (with A c ~ 3.3 in d — 1, see, e.g., 0J) and a suitable order parameter is 
indeed provided by the average density of active sites (n x ) in the stationary state. 

Extensive Monte Carlo simulations of the CP have provided a very detailed quantitative 
description of its dynamical behavior and of the associated non-equilibrium phase 
transition, which turns out to be in the same universality class as directed percolation 
and therefore its universal features are properly captured by the so-called Reggeon field 
theory (see, e.g., E3 E] ) . 

The directed percolation (DP) universality class has been recently the subject of renewed 
interest in the context of ageing phenomena |H1 EH E] > which are characterized by 
two-time quantities (such as response and correlation functions, see below) which are 
homogeneous functions of their time arguments and do not display the time-translational 
invariance which is expected in a stationary state. The two-time quantities of interest 
are the connected density- density correlation function C x ^ x r(t,s) := (n x (t)n x i(s)) — 
(n x (t))(n x >(s)) and the response function R x - x /(t, s) : = S(n x (t)) /8k x >(s)\ k =o where k x >(s) 
is the field conjugate to n x /(s) — corresponding to a local spontaneous activation rate 
of the lattice site x' at time s — and (...) stands for the average over the stochastic 
realization of the process. 

In [HI EI the ageing behavior of the critical CP was investigated numerically. A system 
with a homogenous initial particle density (corresponding to the active phase with 
A ^> A c ) was quenched at t — to the critical point A = A c and R x =o(t, s) and 
C x= o(t, s) were determined, finding the following scaling behavior (t > s) 



R x=0 (t, s) = s- l ~ a f R (t/s), C x=0 (t, s) = s- b f c (t/s), 



(1) 
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where 



f R (y » 1) ~ y 



fc(y > i) ~ y 



(2) 



and b = 2(3/ {vz) = 25 [(3, v, z and 5 := (3/{yz) are the standard critical exponents of 
DP, see, e.g., jl] - the values of z and 5 for d — 1,2 and 3 are reported in table 
In particular it was found that Xc/z = 1.85(10) and Xr/z = 1.85(10) for d — 1 [1113, 
whereas Xc/z = 2.75(10) and Xr/z = 2.8(3) for d = 2 Jf)]. These results suggest that 
Xr = Xc independently of the dimensionality d = 1 and 2. On the same numerical 
footing it was noticed — with some surprise — that 1 + a = b, in stark contrast to the 
case of slow dynamics of systems with detailed balance such as Ising ferromagnets, for 
which a = b (see, e.g., [121 EI]). The fluctuation-dissipation ratio [TH IT2] 



and its long-time limit Xoo := lim^oo lim^oo X(t, s) have been used, for these magnetic 
systems evolving in contact with a thermal bath of temperature T^, to detect whether the 
(equilibrium) stationary state has been reached, in which case the fluctuation-dissipation 
theorem yields = 1. Given that b = a + 1 for the CP, the fluctuation-dissipation 
ratio — as defined in equation © — would always yield a trivial value of Xoo and 
therefore it would not serve its purpose. In [Sj, it was suggested to consider, instead of 



and define Soo := lim s ^oo lim^oo 2(t, s), which has now a finite non-trivial value 
= 1.15(5) for the CP in d — 1 (We comment on the meaning of 5 in 

subsection 12.21 ) In addition, the numerical results of [HI support the applicability 
of the so-called theory of local scale invariance (LSI) [13 E] to the ageing behavior 
of the CP. This theory tries to use local space-time symmetries to constrain the form 
of /i?(t/s) and has been applied to ageing phenomena in magnetic systems quenched 
from an initial high-temperature state to and below the critical temperature pHl E| and 
bosonic reaction-diffusion systems [IB]. Recent simulations suggest that LSI holds also 
in the parity conserving universality class jTHj. However, in the case of critical ageing 
of the Ising model with purely relaxational dynamics, it was shown that corrections to 
the prediction of LSI are present both in field-theoretical results at two loops [T3J [2D] 
and in simulations [2U E2]- A more general version of LSI fT] - the version we shall 
refer to in this paper — improved considerably the agreement with simulations while 
the disagreement with field-theoretical predictions remained. Also in the case of the 



X(t, s) := T b R x=0 (t, s)/d s C x=0 (t, s) 



(3) 



X, 



E(t,s): 



R*=o(t,s) = f R (t/s) 
Cx=o(*,s) fcit/s) 



(4) 
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CP, numerical results indicate that the scaling function predicted by LSI is incorrect 
for t ~ s [TUl E]. It was suggested in [TTJ |2H| that one could possibly account for 
this discrepancy by extending LSI to include also the case of a nonvanishing average of 
the order parameter (in the present version a vanishing order parameter is implicitly 
assumed). As this extension is at present still lacking, in what follows we shall refer to 
the available version of LSI. 

The previously mentioned (numerical) works leave some questions open: 

(i) Could the relations 1 + a = b and \r = Ac have been expected? 

(ii) Is Xr (= Ac) an independent critical exponent? 

(Hi) Is a universal quantity like in ageing systems with detailed balance? 

(iv) To what extent does LSI describe some of the features of the ageing behavior in 



The aim of this note is to answer these questions by adopting a field-theoretical approach. 

The rest of the paper is organized as follows. In section El we introduce the field- 
theoretical model (Reggeon field theory), the formalism and we set up the general 
framework for the perturbative expansion. In addition we discuss the expected scaling 
forms for the two-time response and correlation functions, the relation among the 
different ageing exponents and between X and 5, providing complete answers to the 
questions (i), (ii) and (Hi). In sectional we calculate the response function and its long- 
time ratio Sqo to the correlation function up to first order in the e-expansion (e = 4 — d) 
and then we compare our results to the predictions of LSI, providing an answer to 
question (iv). In section 0] we summarize our findings and present our conclusions. 
Some details of the calculation are reported in the appendix. 

2. The field-theoretical approach 

As explained in detail in [21 IS], the universal scaling properties of the CP (more generally, 
of the DP universality class) in the stationary state are captured by a Reggeon field- 
theoretical action S. In the critical case it reads 



where the <p(x, t) and h(x, t) are the coarse-grained versions of the particle density n x (t) 
and of the spontaneous activation rate of lattice sites n x (t) (external perturbation), 
respectively, whereas <f(x,t) is the response field, u > the bare coupling constant of 
the theory and D the diffusion coefficient (D = 1 in what follows, unless differently 
stated). 



the CP? 
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In terms of the action (|5|). the average of an observable & over the possible stochastic 
realizations of the process is given by {&) = J[dLpdip] 6 q~ s \v>,01^ As a result, the response 
of {&) to the external perturbation h can be computed as 5(0)/5h(x,t) = (<p(x, t)&), 
leading to the following expression for the (linear) response function R x _ x r(t,s) : = 
5(<p(x,t))/5h(x', s)\ h=Q = ((p(x,t)<f(x', s)) h=0 . 

The action fl^J) with h = and t G (— oo, oo) is invariant under the duality transformation 
<f(x, t) * — > -<f(x, -t) (6) 

(the so-called rapidity reversal — RR) which implies that the scaling dimensions [. . .] sca i 
of the fields (p and (p are equal (see, e.g, 0IH]) : 

[cp(x, t)] sca i = [<p(x, t)] sca i = p/v . (7) 

Note that, as in the case of systems with detailed balance (DB), the scaling dimensions 
of the fields (with t > 0) do not change if RR (alternatively, DB) is broken by the 
presence of initial conditions (at time t — 0). 
If RR is a symmetry of S, then 

C x - X >(t,s) := (<f(x,t)(p(x',s)) h=0 - ((f{x,t)) h=0 (<f(x',s)) h=0 

^= (<p( x i -t)<p( x 'i s))h=o - (<p( x i -t))h=o(<p(%', s))h=o = , (8) 

where the last equality is a consequence of causality |21j, and therefore correlations 
vanish in the stationary state of the CP. On the other hand, a suitable initial condition 
(say, at time t — 0) can effectively generate correlations which decay for t > 0. Generally 
speaking one expects this decay to be exponential in t for A ^ A c , due to a finite 
relaxation time both in the absorbing and active phase, whereas an algebraic decay 
is expected in the critical case A = A c we are interested in. (This picture is indeed 
confirmed by the numerical results of [8 .) Accordingly, C x=0 (t, s) := (f(x, t)ip(x, s)) — 
(<f(x, t)) (<f(x, s)) does no longer vanish during the relaxation from the initial condition 
and for t, s ^ it has the same scaling dimension as R x= o(t,s). Two straightforward 
consequences of this fact are: 

(A) 1 + a = b [see equation (JTJ], which could have been expected on this basis and 
therefore is valid beyond the mere numerical coincidencet and 

(B) E(t, s) [see equation - and therefore — is a ratio of two quantities with 
the same engineering and scaling dimensions and therefore it is a universal function 
(see, e.g., j2Ej) which takes the same value in all the models belonging to the same 
universality class (in particular, the lattice CP, the DP and the model (jSJ) on the 
continuum). This is analogous to the case of in systems with DB. 

f The same conclusion can be drawn from the analysis presented in j2H] (see also 
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2.1. Scaling forms 

The non-equilibrium dynamics of the CP after a sudden quench (at time t = 0) 
from the active phase A ^> A c (characterized by {n x (t < 0)) = 1 on the lattice and 
(<p(x,t < 0)) / on the continuum) to the critical point A = A c is partly analogous 
to the non-equilibrium relaxation of Ising systems with dissipative dynamics after a 
quench from a magnetized state (e.g., a low-temperature one with T <C T c , where 
T c is the critical temperature) to the critical point T = T c [27]. In both cases the 
order parameter m(t) := ((p(x,t))$ § provides a background for the fluctuations, with a 
universal scaling behavior (2HI EH] 



where m := (f{x,t = 0)) is the initial value of the order parameter, a is related to the 
scaling dimensions of the fields if and (p in real space via a = (d — [(p] BC a\ — [<p]seai)/z 
(we recall that [y?] sca i = P/v), 9 is the so-called initial-slip exponent and q = 
9 + a + f3/(uz). In equation Tm{v) is a universal scaling function once the non- 
universal amplitudes A m and B m have been fixed by suitable normalization conditions, 
e.g., ^m(O) = 1 and Tm{v — * oo) = v^ 1 + 0(v~ 2 ). [Note that Fm{v oo) ~ v^ 1 
is required in order to recover the well-know long-time decay of the order parameter 
m(t) ~ Within the Ising universality class a — (2 — n — z)/z and 9 P§1IT3] 

whereas within the DP class a = —rj/z (indeed [v?] sca i = fiscal = (d + rj)/2 0) 
and 9 = |2H]- In both cases it turns out that the width A = ([(p(x,0) — m } 2 ) 
of the initial distribution of the order parameter (p - - assumed with short-ranged 
spatial correlations — is irrelevant EH] and controls only corrections to the leading 
scaling behavior, so that an initial state with A 7^ is asymptotically equivalent to 
one with no fluctuations Ao = 0. Equation (j^J clearly shows that a non- vanishing 
value of the initial order parameter introduces an additional time scale in the problem 
r m = (B m mo) -1 ^ P3I2H1I2I!- The long-time limit we are interested in is characterized 
by times 3> r m and therefore the relevant scaling properties can be formally explored in 
the limit m — > 00, i.e., r m — ► 0. 

According to the analogy discussed above one expects the following scaling form for the 
two-time response function in momentum space (see, e.g., EH] and section 3 in fI7\ ) 



I For magnetic systems tp represents the local magnetization whereas in the present case (tp) is the 
coarse-grained density of active sites. 

§ Hereafter we assume invariance of the system under space translations. 



m(t) = A m m Q t e+ ^ M (B m m t') 



(9) 



R q=0 {t >s,s)= A R (t - s)\t/s) e F R (s/t, B m ml'H) , 



(10) 
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where Fr is a universal scaling function once the non-universal amplitude Ar has been 
fixed by requiring, e.g., F R (0,0) = 1. In the case of the DP universality class, the 
two-time Gaussian correlation function Cx(t,s) vanishes for itlq = — apart from 
irrelevant corrections due to a finite A — whereas the two-time Gaussian response 
function R x °\t, s) is invariant under time translations. Taking into account causality, it 
is easy to realize that all the diagrammatic contributions to the full response function 
(with M ^ 0) have the same invariance and therefore F R (s/t,0) = 1. (As opposed to 
the case of the Ising universality class, where 9 ^ and Fu(s/t,0) ^ 1 (SOI EH EH EI]-) 
In fact, the case m = would correspond on the lattice to an initially empty system 
with A = A c where nothing happens as long as one does not switch on the external field 
K x /(s). Accordingly, the response to k x /(s), at time s + At has to depend only on At. 

In the long-time limit t > s ^> r m one expects Fr to become independent of the actual 
value of m and therefore 

F R (x,v^oo) = ^a*- 9 F R (x), (11) 

where 9 is an additional exponent [related to Xr, cf equation p8jl ]. so that the resulting 
scaling is 

R g=0 (t >s,s)= A R (t - sf{t/s) S T R {s/t) . (12) 

Ar is a non-universal constant which can be fixed by requiring ^^(0) = 1 and which 
has a universal ratio to Ar [see equation For the Ising universality class, is was 

found gZj 

$ = -(l + a + ^j. (13) 

On the other hand, it is easy to see that the scaling arguments which were used to drawn 
this conclusion for the response function (see section 3 in |2Zj) are valid also for the CP. 

In the scaling behavior of the two-time correlation function C q= o(t, s), after a quench 
from the active phase to the critical point, has been discussed with the result that (see 
equations (3.3) and (4.14) in [25J) 



C q=0 {t >s,s) = A c {t - sf{t/s) e Tc{s/t) (14) 

where the non-universal amplitude Ac can be fixed by requiring J^c(0) = 1, and 9 is 
indeed given by equation ()13j). Note that in analogy with S [see equation (J1J)] one can 
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g, v ._ gg=ofa s) _ A R F R {s/t) 
J ' C q=0 (t,s) AcTc(s/t) 

[where we have used the scaling forms ()12|) and ()14|) ] and 

lim lim H(i, s) = — — (16) 



s— >oo t— >oo ' ^4, 



c 



which S(i, s) and Soo, a universal function and amplitude ratio, respectively. 

Although in general E(t, s) ^ E(t, s), the argument discussed in [3UJ HSj can be used to 



conclude that = also in this case. 



In what follows we focus on R q= o and C q= o, i.e., the response and correlation function of 
the spatial average of the density of active sites (N~ 1 (^2 X n x ) on a lattice with N sites). 
On the other hand, the corresponding scaling forms equation (I12j) and (fl"4"j) can be easily 
generalized to q ^ taking into account that this amounts to the introduction of an 
additional scaling variable y = ArjDq z (t — s) where Ad is a dimensional non-universal 
constant which can be fixed via a suitable normalization condition. Accordingly, the 
scaling behavior of the autoresponse R x =o and autocorrelation C x= o functions can be 
easily worked out from equations (j!2)l and (|14)1. leading to the identification of the 
exponents a and A^c in equations (JTJ) and © as 

l + a = b = --a=^- = 25 (17) 

z vz 

and 

— = — = a — = l + d-\ — . (18) 

z z z z 

Therefore we conclude that: 

(C) Xr is equal to A^ beyond the mere numerical coincidence observed in jHUH] and 

(D) Xr = Ac is not an independent critical exponent but it is given by 

X R = z + zd + d. (19) 

In table ^ we compare the values of Xr^c obtained from this scaling relation and 
the available estimates of 5 and z to the results of fitting the asymptotic behavior 
of the scaling functions fc,R [see equations (JTJ) and (J2J)] from numerical data [HUH]- 
The results reported in the last three columns are in quite good agreement both 
for d — 1 and 2. 
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a 


c 




z 


1 + o + d/z 


Xr/z 


Xc/z 










1.85(10) |8J 


1.85(10) jHT 


1 


0.159464(6) 


1.580745(10) 


1.792077(10) 


1.9(1) i 

1.76(5) noi 


1.9(1) i 


2 


0.451 


1.76(3) 


2.58(2) 


2.75(10) [9 


2.8(3) P 


3 


0.73 


1.90(1) 


3.30(1) 


not available 



Table 1. Comparison between the direct numerical estimates of \r,c / z and the 
predictions of the scaling relation X R /z(— \c/z) = 1 + 5 + d/z [see equation ifTSjl] in 
which the available estimates of S and z (taken from table 2 in 0]) are used. The values 
of ^R,c I z reported in the last two columns have been obtained by fitting the scaling 
behavior of the autoresponse and autocorrelation function determined via density- 
matrix renormalization-group computation [Hj and by Monte Carlo simulations O I10| 
of the contact process. 

The conclusions (A-D) we have drawn so far provide complete answers to the questions 
(i-iii) which we have posed at the end of the Introduction. 

2.2. Relation between X and H 

In [B] H has been introduced on the basis of a formal analogy with the fluctuation- 
dissipation ratio X(t,s). Here we argue that indeed X and 5 play the same role in 
different circumstances. 

The stationary state of a system with detailed balance is characterized by the time- 
reversal (TR) symmetry of the dynamics, which qualifies the state as an equilibrium one 
with a certain temperature T. Time-translational invariance (TT) and TR symmetry 
of the dynamics in the equilibrium state imply the fluctuation-dissipation theorem 
TR x>q (t,s) = d s C x>q (t, s), leading to X = 1 [see equation ©]• In this sense X ^ 1 
is a signature of the fact that the system has not reached its stationary (equilibrium) 
state. (This is usually due to slowly relaxing modes which prevent the system from 
"forgetting" the initial conditions of the dynamics jHHEIj-) When detailed balance does 
not hold, TR is no longer a symmetry of the stationary state. In the specific case of 
the DP universality class, however, the non- equilibrium stationary state is characterized 
by a different symmetry, the rapidity-reversal [RR, see equation (0)], which leads — as 
discussed in subsection 12.11 equation (jHJ) — to vanishing correlations and therefore to 
S -1 , = 0. Accordingly, H -1 , ^ provides a signature of the fact that the system 
is not in its stationary state. In this sense S (S) is analogous to X for both of them 
indicate if the stationary state is eventually reached. In particular this occurs generically 
after a perturbation (e.g., a sudden change in the temperature T or the spreading rate 
A) and therefore X^ = 1 or E^ 1 = = 0. On the other hand, it may happen that for 
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particular choices of the external parameters (e.g., T or A), slow modes emerge which 
prevent the system from achieving its stationary state. As a consequence, during this 
neverending relaxation of the system, TT is broken together with the possible additional 
symmetries which characterize the stationary state. || This is what happens when ageing 
takes place in critical systems with detailed balance. In the next subsection we explicitly 
show that this is also the case for the contact process. 



2.3. The Gaussian approximation 



The analytic computation of the response function follows the same steps as those 
discussed in [2Z] (see also for the non-equilibrium relaxation of the Ising model 
relaxing from a state with non-vanishing value of the magnetization. Here we briefly 
outline the calculation. We assume that the initial state has a non-vanishing order 
parameter mo = (<f{x,t = 0)) which is spatially homogeneous, so that the ensuing 
relaxation is characterized by translational invariance in space. It is convenient to 
subtract from the order parameter field (fi(x,t) its average value m{t) = (<p(x,t)), by 
introducing 

ijj{x, t) :— <p(x, t) — m(t) and ip(x, t) := (p(x, t) , (20) 

leading to (ip(x,t)) = during the relaxation. The action 5* [see equation (jSJ)] becomes, 
in terms of these fields, 



S[if>, #] = I d d xdt |$ [d t - V 2 + 2a] ip - a^ 2 - u(tp - , (21) 

where a{t) := um{t) satisfies the equation of motion (no-tadpole condition) 

d t a + a 2 + u 2 (i) 2 ) = 0. (22) 

Equations and (}2*T]) are the starting point for our calculations. The Gaussian 

response and correlation function are obtained by neglecting anharmonic terms in 
equation (}2*T]) : 

Rf\t,s) = 0{t- s)expl^-q 2 {t- s) -2 dtV(t')j and (23) 

Cf(t >S ) = 2 ^dt'RfMo-^Rf^t') (24) 



o 



|| The one described here is the typical pattern of the spontaneous breaking of a symmetry. See, 

e.g., m 
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where t < := min{t,s}. [In what follows we denote the results of the Gaussian 
approximation with the superscript (0).] In particular, solving (}2~2~)) with u = one 
finds (see also [2H] and the appendix of ~fl|) a^°'(t) = (t + o-q 1 )^ 1 , which scales according 



to © (/?(°) = 1, 
Bm = u and T^ya 



,(0) 



1/2, 2 



(0) 



2 and 77 



(0) 



0, see, e.g., [24J) with A 



(0) 



t > s > r„ 



(B m m ) 



[1 + x) As explained in section l2~T| the leading behavior for 



(J 1 can be explored by taking the limit m oc a — > 00 



from the very beginning. Accordingly, cr^(t) 

2 



i 1 and 



4 0) (t, 



0(t-a 

\ i / 



<z 2 (*-s) 



diW 2 *', 



(25) 
(26) 



which (for q = 0) display the expected scaling behaviors (fT2~j) and (jH|) with ^4 



(0) 
c 



1/2 and .F ( ' 



(o) 

C,R 



X 



1. Equations (USD and (USD yield H^(t, s 



(0) 

7(0) 



which is exact for d > 4 but it is almost twice as large as the result which was found 
in d — 1, = 1.15(5) jH]. As explained in subsection 12.11 7^ is a signal of the 
breaking of the RR symmetry characterizing the non-equilibrium stationary state of the 
CP (and of the DP universality class in general). Here we show explicitly (within the 
Gaussian approximation — though the conclusion is expected to be valid also beyond 
the approximation) that only the homogeneous fluctuation mode q = is, at criticality, 
responsible for such a breaking, as in the case of systems with detailed balance (see, 
e -g> [EH 123 EH] ) • Indeed, let us generalize equation (|T5j) to modes with q 7^ 0: 



^(t > s,s) :-- 



C q (t,s) 
R q (t,s) 



(27) 



Within the Gaussian approximation one readily finds [see equations ()25|) and ([26p] 



?(o) 



(t,s) 



14! 
2^ 



-1+2/- 



y_,y_ 

2! 3! 



(2f 



y=2Dq 2 s 



where the additional scaling variable y = A£,Dq z s = Dq 2 s appears. This expression 
- as the Gaussian fluctuation-dissipation ratio for systems with detailed balance (see, 
e.g, [30 ) - - is independent of t and depends on y only. In particular, in the long- 



time limit s ^> q one has 



S(o) 



(y » 1; 



-1 



l/y — > 0, for any mode with q 7^ 0, 



indicating that the RR symmetry is asymptotically realized. On the other hand, for 



q 



0. 



1/2, independently of s. For a quench into the active phase 



A > A c the action S and the propagators (|23j) and (|2*S)l get modified according to 
q 2 — * q 2 + r, with r > (see, e.g., 121 )■ Therefore in this case, y = 2D(q 2 + r)s, 



yielding 



(y » 1) 



for s ^> r 1 and independently of q, in agreement with 



what has been stated in subsection 12.21 and with the available numerical evidences jHj. 
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0.5 - 



CO 

o 

<[I] 



0.4 A 



0.3 



0.2 



0.1 



1 1 1 1 1 1 1 ■ 1 1 1 

10 20 30 40 50 

y = 2Dq 2 s 

Figure 1. Ratio [Eg (t, s)]" 1 of the correlation function to the response function [see 
equation 1|27[) ] within the Gaussian approximation (which becomes exact for d > 4). 
In the long-time limit s — > oo, (t, s)] _1 vanishes, unless q = 0, in which case it 
takes the value 1/2. [S^] -1 ^ signals the breaking of RR [see equation © and 
subsection \2.'2\ . 



3. The response function 



In this subsection we determine the response function up to first order in e = 4 — d (4 
being the upper critical dimension d u d of the model [Sj, above which the Gaussian results 
become exact). Some of the details of the calculation are provided in the appendix. 

For future reference we recall that the critical exponents of the DP universality class 
are 0SE1 rj = -e/6 + 0(e 2 ), z = 2 - e/12 + 0(e 2 ), fi/v = {d + r))/2 and therefore [see 
equation (fTTTf> ] 



. r L = j_j r Qi^\ —a a— O i 6 /"I/. -2 



12 



and 6 = -2 + - + 0(e 2 ). 
6 



(29) 



The expression for the renormalized response function is (where x := s/t < V 



RR, q =o(t, s) 



.i- < 1 r 

.2 k.m , 



1 , 1 , 

— m x H In s 

4 12 



7T ln(l — x) 1 T . . . 11 , , x 

1 U 2 (x) n(l - x) H 

12 x 2 v ; 12 v ; 12 



+ 0(6 2 ), (30) 



(see the appendix for details) which can be cast in the expected scaling form (|12p with 
A R = l-e[l-'Q+0{ t ?), (31) 
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F R (x) = l + e 1 + -: + 



1 



l] ln(l - x) - ^Li 2 (x) + 0(e 2 ), (32) 



x 



[note that J^(0) = 1 and T R {x -> 1~) = 1 + e(13 - tt 2 )/12 + 0(e 2 )] and the proper 
exponents (|29p. 

5. i. Comparison with LSI 

For the response function, LSI provides the following prediction jllj : 



where x : = s/t and A, B, C and ro are free parameters. In a first version of LSI it was 
assumed A = C but more recently it has been suggested that this constraint can be 
relaxed [HIE]- On the other hand, requiring R q= o to have the correct scaling dimension 
(see subsection 12. lj) . leads to A = a = —25 + d/z. In addition, the well-established short- 
time behavior s/t — > 0, from equation (j!2|) . leads to B = —9 — a = 1 + 5. The last free 
parameter, C, can be in principle fixed by requiring the response function to display 
a quasi-stationary regime (analogous to the quasi-equilibrium regime in systems with 
detailed balance, see, e.g., Section 2.5 in ^Sj) for At <^ s where At := t — s, in which the 
behavior of the system depends only on At and not on s. This would require C = A [10J. 
In turn, comparing with equation (|T2"j) . this would imply Fr(x) = 1, Vx G [0, 1] which is 
in contrast with the analytical expression in equation (|3*2"|) . The discrepancy between the 
actual scaling behavior of the CP and the prediction of LSI with A = C, was originally 
overlooked in [H] but then it became apparent for t/s < 2 in the detailed numerical 
analysis carried out in [THj, which is also supported by our analytical result. On the 
other hand, if one questions the existence of such a quasi- stationary regime in the CP, C 
remains an available fitting parameter. By adjusting it the agreement between numerical 
data and LSI improves significantly JT], down to t/s — 1 ~ 10" 1 , being discrepancies 
observed only for smaller t/s — 1 [TQJ . It has been suggested [TT] that such discrepancies 
might be caused by corrections to scaling due to having t — s of the order of microscopic 
time scales when exploring the regime t/s — 1 ~ with finite s — However, this does 
not seem to be the case for the data reported in [TU] . 

In order to compare the analytic expression of R q=0 (t,s) to the prediction of LSI 
with C 7^ A we introduce the quantity TZ(x) := A]^(t — s)~ a (t/s)~ d R q= o(t, s) where 
Ar is fixed by the requirement TZ(x —>■ 0) = 1. Accordingly, the field-theoretical 
prediction (|T2"|). yields (for d < 4) 



R q =o(t > s, s) = r s A x B (l - x) c , 



(33) 



^ft(x) = F R {x) = 1 + e f(x) + 0(e 2 ) 



(34) 
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where f(x) is given by the expression in square brackets in equation (|52Jl : 



f(x) := 1 + 



x 
12 



(H 



ln(l -x) - -U 2 (x), 



(35) 



with f(l) = (13 - tt 2 )/12. LSI predicts, instead, 



^lsi(^) = (1 - ar" 



AC 



(36) 



where AC := C — A. Above the upper critical dimension d uc d = 4, 11ft (x) = 1 and 72-lsi 
fit it with AC = 0. Assuming continuity of critical exponents, one expects for d < 4 (i.e., 
e > 0), AC = c\€ + 0(e 2 ), — where c\ is a fitting coefficient — and therefore T^lsi(^) — 
l + eciln(l-x) + 0(e 2 ). Of course ft FT -ft LS i = e[f(x)-c a ln(l-x)] + 0(e) ^ for e > 
(i.e., d < 4) whatever the choice of c\ is, as it is clear by inspecting the expression (}35|) . 
Strictly speaking LSI in its present form is not a symmetry of R q= o{t, s), whatever AC 
is, at least sufficiently close to d nc d (small e). Nonetheless, the shape of f(x) is such that 
a proper choice of c\ can reduce the difference H-ft — in figure El we report the 

comparison between Xvl'Rft^x) = tf( x ) + 0{e 2 ) and lnT^LSi^) — ec i l n (l — x) + 0(e 2 ) 
with ci = (i.e., C = A) and ci = —0.081(2) (corresponding to C ^ A). This plot 



«3 




Figure 2. Comparison between e 1 ln7?.FT(^) and e 1 In 7?.lsi (x) for two different 
values of the fit-parameter c\ . 



clearly indicates that the scaling function 1Zft{%) shows corrections to the behavior 
predicted by LSI with C = A already for t/s — 1 < 8, i.e., 0.11 < x < 1, for which 
7e FT (x)/^ LS i(x) > l+0.01e+O(e 2 ). On the other hand the same discrepancy is observed 
only for t/s-1 < 0.2[0.13], i.e., 0.83[0.88] < x < 1 for Cl = -0.083[-0.079] and therefore 
a suitable choice of c\ (i.e., C) reduces by more than one order of magnitude the value 
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of t/s — 1 below which the difference between the predicted scaling behavior and LSI 
becomes sizable. This is analogous to what has been observed in numerical data fUJ E2] 
for the scaling of the autoresponse function R x=Q (t,s). (See, e.g., figure 2 in [10].) 
On the other hand, the field-theoretical results presented here are free from numerical 
artefacts and corrections to scaling, which were invoked in ^T] to explain the findings 



As a result of our fitting procedure, we have determined the "e-expansion" of the 
exponent AC, finding AC = — 0.081(2)e + 0(e 2 ). According to the notations adopted 
in [TT] AC = —(a' — a) and therefore we conclude that a' — a = +0.081(2)e + 0(e 2 ). 
For the one-dimensional contact process it was found a' — a = +0.270(10) ^IJ O-O %, 
whose sign is indeed in agreement with our result (generally speaking, the sign of the 
correction to the Gaussian behavior is correctly captured by the lowest non-trivial order 
in the standard e-expansion). Insisting in interpreting our result as an e-expansion we 
can even attempt an estimate of a' — a in physical dimensions e = 1, 2 and 3. It comes as 
a surprise that not only the order of magnitude is the correct one but also the estimate 
for the one-dimensional CP (i.e., e = 3) a' — a ~ +0.24 is in very good agreement with 
the actual numerical result. 

In passing, we note that a mechanism similar to what has been described here could 
possibly explain the (apparent) agreement between LSI with C ^ A and the scaling form 
of the integrated response function of the magnetization in the Ising model quenched 
from high temperatures to the critical point [2*T| 111]. 

3.2. One-loop prediction for H 

The correlation function was already calculated in [25], up to one loop in the e-expansion, 
with the result 



of no]. 



2 A 



Ax 1 



€ F(oo)]+0(e 2 ) = l + 



e /9vr 2 



361 



) 



+ 0(e 2 ) (37) 



c — 



6 V 20 



80 



and 




(38) 



% Note that in JU| the value of a' 



is reported with the wrong sign. 
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where A and F are given, respectively, in equation (4.4) and (3.5) therein. Accordingly, 
we obtain for the ratio in (jl6j) the value 

+ 0(e 2 ). (39) 

Clearly the Gaussian value of gets a downward correction (119/480— 7r 2 /120 ~ 0.166) 
for e > 0. This is qualitatively satisfying, as it brings the theoretical estimate closer to 
the value which was determined numerically for d — 1 [2^ = 1.15(5)]. In addition, a 
direct extrapolation of equation (f3TJ|) to e = 3 (neglecting possible 0(e 2 ) terms) yields 
S^ 1 = 1.01 which is — for a one-loop calculation — in a surprisingly good agreement 
with the actual value. + More reliable theoretical estimates can only be based on the 
knowledge of higher-loop contributions. 

4. Conclusions 

In this note we have studied the ageing behavior in the contact process by adopting a 
field-theoretical approach. The results presented here equally apply to all the models 
belonging to the same universality class as the contact process, such as the directed 
percolation. 

We confirmed analytically the scaling behavior 

R x=0 (t, s) = s- l - a f R {t/s), C x=0 (t, s) = s- b f c (t/s), (40) 

which was observed numerically for the two-time response and correlation functions [HI El 
after a quench from the active phase to the critical point. In addition we have shown 
that: 

(a) The relation b = 1 + a, due to the rapidity-reversal symmetry, and Xr = Ac [see 
equation (J2J)], hold beyond the numerical evidences provided in [HUH]- 

(b) \r,c are related to known static and dynamic exponents via Aj^c = z(l + 8) + d [see 
equation ()19|)]. which is confirmed by the available numerical data (see table Q). 

(c) The ratio 5 of the response to the correlation function [equation (j3J)] defines a 
universal function and its long-time limit a universal amplitude ratio, which 
can therefore be studied within a field-theoretical approach. 

The result of the e-expansion around d = 4 [see equation ()39|)] is in qualitative 
agreement with the numerical estimate of S^ 1 [H] . 

+ Extrapolating the result (S^ as S^ xtr) = 2/[l + e(119/480 - n 2 /120)} would give 1.33 for e = 3, still 
not too far from the numerical result. 



119 
480 



7T" 

120 
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(d) After a quench, S^ 1 ^ for the CP, like — 1 ^ for systems satisfying detailed 
balance, is a signature of the spontaneous breaking of the symmetry (5 associated 
with the corresponding non-equilibrium (6 = rapidity-reversal) and equilibrium 
(& = time-reversal) stationary states. 

(e) The universal scaling function Fr^x) - - describing the ageing behavior of the 
response function R q =o(t, s) [see equation lJI2jl ] — does not agree with the prediction 
of the theory of local scale invariance when the non-Gaussian fluctuations are taken 
into account [specifically already at 0(e) for e > 0, see equation (J52*|) ]. This provides 
an additional example against the general applicability of the present form of LSI 
to critical systems. It remains to be seen whether an extension of LSI to systems 
with nonvanishing order parameter will be able to describe correctly the scaling 
function of the CP. 

Acknowledgments: We thank Tilman Enss, Malte Henkel, Michel Pleimling, Beate 
Schmittmann and Uwe Tauber for useful discussions and comments and the organizers 
of the school N on- Equilibrium Dynamics of Interacting Particle Systems at the Isaac 
Newton Institute of Mathematical Sciences, Cambridge, where this work was initiated. 
FB acknowledges the support by the Deutsche Forschungsgemeinschaft through grant 
no. PL 323/2. 

Appendix A. Diagrammatic expansion of the response function 

In this appendix we provide some details of the calculation of the response function. 
First of all we solve equation (|22[) — to first order in the effective coupling constant 
v := (8n)~ d / 2 u 2 / D 2 — where (?/> 2 ) is given by 



For later convenience we introduce the notation (dq) := d d q/(2n) d [and dimensional 
regularisation of the integrals is understood] and = 12r(l — d/2)/T(5 — d/2) = 
— 12/e + 3 + 0(e). Accordingly, a(t) is given by 



where o~(°'(t) = 1/t is the solution within the Gaussian approximation. (As we 




(A.l) 




(A.2) 



are interested in the long-time limit, we have taken o"o 



oo.) As expected, the 
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(h) 



h 



Figure Al. Diagrammatic contributions to the response function up to first order in v. 
Directed and undirected lines represent response l|25(l and correlation (|26(l propagators, 
respectively. 



expression (|A.2|) has a dimensional pole for e — > 0, which is removed once the 
renormalized quantities (with subscript R) are introduced. In particular (following j3j): 



ipR = Z l J 2 ^, = Z^!^, D R = Z D D, and v R = Z^y, e v, (A. 3) 



A/2 



where \i is an arbitrary momentum scale and Z$ = 1— 4v/e+0(v 2 ), Z D = l+2v/e+0(v 2 ) 
and Zfj = 1 — 2Av/e + 0(v 2 ) are suitable renormalization constants. Taking into account 
the definition of a: a (t) = {%n) d ^Dv ll2 m{t) and that m R = Z^ 2 m } one finds for the 
renormalized function o~ R (t) at the fixed-point v* R = e/24 + 0(e 2 ) [3] 



a R (t) = ~ (l + |lnt-9e) + 0(e 2 



(A.4) 



which scales according the long-time limit of equation © with the proper exponents 
[see before equation in section 0] and A m /B m = [1 — 9e + 0(e 2 )]/u*. To compute 
the response function to first order in v, one has to take into account that the Gaussian 
propagator gets modified because of the 0(v)-term in equation (|A.2|) [27 . Using (|A.2|) 
in equation one gets (hereafter t > s and a; := s/t) 



X 



1 + vs 



24 



In x — 18 In x — 6 In x \ + 0(v ,ev 



.(A.5) 



In addition, there are two diagrams contributing to the response function at order v, 
depicted in figure IA"T1 Notice that we can use the O(u )-propagators to compute them 
up to 0(v). We first need the expressions for the following bubbles 



B RC (t 1 ,s 1 ) := J{dq)Rf\t 1 ,s 1 )C^(t 1 ,s 1 ) 

Jo 



m M fi ) (/-, - ,, )- 



(A.6) 
(A.7) 
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for ti < si, whereas for si > t\ both B RC and B RR vanish because of causality. The 
contribution Ii is then given by the integral 



h 



dh / ds 1 R%(t,t 1 )B RC (t 1 ,s 1 )R%(s 1 ,s) 



(A.8) 



which has essentially been calculated in [21]. (Indeed, taking into account that 
B R c(ti,si) here is equal to £j~ i?d(ti,si) of |2Z| ■ - see equation (A. 4) therein — one 
finds that the rhs of equation (jA.8|) is equal to (s/t) 1//2 J 2 (t, s)/3 where J 2 is given by 
equation (A. 11) of 27J.) The result is 



4, 2tt 2 x /l 
-ma; — 1 + xH h 6 m(l — x) I 1 

e 3 \x 



3 In x 



+ ln 2 x - 4Li 2 (x) + 0(e) 
The contribution J 2 is instead given by 

h:= f dt x ^ ds 1 Rf2 (t,t 1 )B RR (t 1 ,s 1 )Rf2 (s 1 ,s) 



(A.9) 



(87r)"^x 2 I dt±ti 2 f ds 1 s 2 1 (t 1 -s 1 ] 

J s J s 



e/2-2 



(A.10) 



We compute first the integral over s 1; which yields after some manipulations and for 
generic e: 

_ ( tl - a) e/2-l a 2 2fti _ 3 )e/2-l tia 



dsi s\{ti — S] 



+ 



e/2 + 1 (e/2-l)(e/2 + l) 

2(ti - s) e/2 ti 



(A.ll) 



( e /2 - l)e/2(e/2 + 1) 
The resulting three integrals over t\ are easy and one finally gets: 



371") ix 2 S t/2 



2 4 In x 9 7r 2 
1 In x H 

e e 3 



+ 2Li 2 (x) +x-ln(l-x) -2 
i \ _ 



(A.12) 



where the relation Li 2 (l — x^ 1 ) = — ln 2 x/2 — 7r 2 /6 + lnxln(l — x) + Li 2 (x) was 
used. Collecting the different contributions we get the final expression for the response 
function: 



R q= o(t, s) = R'}%(t, s) + Avh - 2vl 2 + 0(v 2 , ev) 



(A.13) 



x 



l + vs*' 2 



- -61nx + 27r - 12Li 2 (x) 



-221n(l -x) +2x + 



241n(l -x) 



x 



+ 0(v 2 ,ev). (A.14) 



Its dimensional pole is properly removed by introducing the renormalized quantities 
according to equation (|A.3|) . so that R R , q =o = Z^R q=0 which gives indeed the expression 
in equation flUD for v = v* R = e/24 + 0(e 2 ). 
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